Singular instability of exact stationary solutions of the nonlocal Gross-Pitaevskii 
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In this paper we show numerically that for nonlinear Schrodinger type systems the presence 
of nonlocal perturbations can lead to a beyond-all-orders instability of stable solutions of the lo- 
cal equation. For the specific case of the nonlocal one-dimensional Gross-Pitaevskii equation with 
an external standing light wave potential, we construct exact stationary solutions for an arbitrary 
interaction kernel. As the nonlocal and local equations approach each other (by letting an appropri- 
ate small parameter e — » 0), we compare the dynamics of the respective solutions. By considering 
the time of onset of instability, the singular nature of the inclusion of nonlocality is demonstrated, 
independent of the form of the interaction kernel. 
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In almost all applications where the nonlinear 
Schrodinger (NLS) equation is relevant, it arises as a 
simplified model of a nonlocal description. This is true 
in water waves, plasma physics, nonlinear optics, and 
Bose-Einstein condensates (BECs). For water waves |l]] 
and plasmas 0, the nonlocality occurs due to a Fourier 
transformation of nonlinear equations. In nonlinear op- 
tics H , the nonlocality is due to the spatial dependence of 
the susceptibility tensor. Locality in each of these cases 
is obtained as a quasi-monochromatic approximation of 
the nonlocal model. In plasma physics the nonlocal effect 
is known as Landau damping. The use of mean-field the- 
ory results in nonlocality for BECs J§-^|, which reduces 
to the Gross-Pitaevskii (CP) equation (the NLS equation 
with external potential) when assuming a hard pairwise 
interaction potential. As such, the consideration of non- 
local perturbations of solutions of the NLS equation is 
important. Such perturbations result in the NLS equa- 
tion with the cubic nonlinearity replaced by a nonlocal, 
nonlinear term. The nonlocal, nonlinear term is a con- 
volution of the modulus squared of the solution with an 
interaction kernel, prescribed by the physical problem. 
In this paper our objective is to numerically examine the 
stability of solutions of this nonlocal NLS equation, and 
discuss how it differs from the local description. 

To consider a specific nonlocal model, we examine 
mean-field theory of many-particle quantum mechanics 
with the particular application of BECs trapped in a 
standing light wave. The classical derivation given here 
is included to illustrate how the local and nonlocal mod- 
els are related. The inherent complexity of the dynamics 
of N pairwise interacting particles in quantum mechanics 
often leads to the consideration of such simplified mean- 
field descriptions. These descriptions are a blend of sym- 
metry restrictions on the particle wave function 0] and 
functional form assumptions on the interaction potential 
Here we do not impose any assumptions on the 
pairwise interaction potential. 



The dynamics of N identical pairwise interacting quan- 
tum particles is governed by the time-dependent, A^-body 
Schrodinger equation 



dt 2m 



N 

-E 

<=i 



Wixi-x^ + Y^Vixi)*, (I) 
i 



where = (xu, x i2 ,x i3 ), * = $>(xi,x 2 ,x 3 , ...,x N ,t) 
is the wave function of the Af-particle system, A N = 



J2 

J Xi2 



is the kinetic energy 



or Laplacian operator for A^-particles, W(xi — Xj) is the 
symmetric interaction potential between the i-th and j- 
th particle, and V(xi) is an external potential acting on 
the i-th. particle. Also, Ti is Planck's constant divided by 
2-7T and m is the mass of the particles under consideration. 

One way to arrive at a mean-field description is by 
using the Lagrangian reduction technique pj , which ex- 
ploits the Hamiltonian structure of Eq. ([l]). The La- 
grangian of Eq. ([!]) is given by Q 
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-J2(W(x l -x J ) + V(x l ))\y\ 2 \dx 1 ---dx N (2) 
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The Hartree-Fock approximation (as used in ji) for 
bosonic particles uses the separated wave function ansatz 



* = ■tp 1 (x 1 ,t)ip 2 (x 2 ,t) ■ ■■ip N (x N ,t) 



(3) 



where each one-particle wave function ip(xi) is assumed 
to be normalized so that (ip(xi)\ip(xi)) 2 = 1. Since iden- 
tical particles are being considered, 



(4) 



enforcing total symmetry of the wavefunction. Note that 
for the case of BECs, assumption (|J) is only approximate 
if the temperature is not identically zero. 
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Integrating Eq. (g) using (g) and (H) and taking the 
variational derivative with respect to ipfei) results in the 
Euler-Lagrange equation 



Hi- 



dip(x, t) 



2m 



Aip(x,t) + V(x)ip(x,t) 



+ (N-l)ip(x,t) W(x-y)\Tjj(y,t)\ 2 dy. (5) 

J — oo 

Here, x = Xi, and A is the one-particle Laplacian in 
three dimensions. The Euler-Lagrange equation (||) is 
identical for all ip(xi,t). Equation (||) describes the non- 
linear, nonlocal, mean-field dynamics of the wave func- 
tion ij){x,t) under the standard assumptions (||) and (^) 
of Hartree-Fock theory The coefficient of ip(x,t) in 
the last term in Eq. (||) represents the effective poten- 
tial acting on ip(x,t) due to the presence of the other 
particles. 

At this point, it is common to make an assumption on 
the functional form of the interaction potential W(x — y). 
This is done to render Eq. ([5]) analytically and numeri- 
cally tractable. Although the qualitative features of this 
functional form may be available, for instance from ex- 
periment, its quantitative details are rarely known. One 
convenient assumption in the case of short-range poten- 
tial interactions is W(x — y) = kS(x — y) where 6 is the 
Dirac delta function. This leads to the Gross-Pitaevskii 
HU mean-field description: 
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where f3 = (N — 1)k reflects whether the interaction is 
repulsive (/3 > 0) or attractive (/3 < 0). This assump- 
tion on the interaction potential W(x — y) is difficult to 
physically justify. Nevertheless, Licb and Seiringer ||] 
show that Eq. (0) is the correct asymptotic description 
in the dilute-gas limit. In this limit, Eqs. and (£) 
are asymptotically equivalent. Thus the nonlocal Eq. (6) 
can be interpreted as a perturbation to the local Eq. (j^). 
Note that the results of || do not have implications for 
the asymptotic equivalence of the stability of solutions. 

Since their first successful demonstrations in 
1995 continuous progress is being made in trap- 

ping, controlling, and manipulating Bose-Einstein con- 
densates in a variety of experimental configurations . 
Although many experiments rely solely on harmonic 
confinement to trap the condensate, we consider the 
situation of an external standing-light wave potential 
within a confining potential |l2|-|l4|] . This standing- light 
wave pattern is generated by the interference of two 
quasi-monochromatic lasers in a quasi-one-dimensional 
configuration. The quasi-one-dimensional regime holds 
when the transverse dimensions of the condensate are 
on the order of its healing length and its longitudinal 
dimension is much longer than the transverse ones. The 



rescaled governing mean-field evolution (j^) in the quasi- 
one-dimensional regime is given by 

dib 1 d 2 ib f°° 
i^ = --Q^+ai>J R(x - y)my 7 t)\ 2 dy + V(x)^. (7) 

Here a = ±1 is the sign of the interaction potential 
W(x — y) at close range. Thus, a determines whether 
the close-range interaction is repulsive (a — 1), or at- 
tractive (a = —1). Hence, a — sign(a), where a is the 
s-wave scattering length of the atomic species. Depend- 
ing on the species, a is either positive or negative, so 
that both signs of a = sign(a) are relevant for BEC ap- 
plications. With these defintions, R(x — y) is the rescaled 
interaction potential, which is positive at close range and 
is normalized to unity J_ R(z)dz = 1. This normaliza- 
tion condition is equivalent to a rescaling of variables. 

The external potential which models the standing light 
wave is given by 0,H 



V(x) = V sin 2 (fca;) 



(8) 



where k is the wavelength of the periodic potential. 

The nonlocal, nonlinear equation (Q) with the periodic 
potential (||) admits a one-parameter family of exact so- 
lutions. These solutions are found using an amplitude- 
phase decomposition 

t/j(x, t) — r(x) exp [i9(x) — iwt] . (9) 

Then 

r(x) 2 = Asm 2 (kx) + B 
where B is a free parameter and 
-V 
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A(k) = 
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tan(#(x)) = Wl - 
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tan(fcc), (lib) 



Lo(k) = z 1- aB 
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/OO 
R(z) cos(2fcz) dz = R(2k), (lid) 
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where R(2k) is the Fourier transform of R(z) evaluated 
at 2k. Equations ( pi] ) can be verified by direct substitu- 
tion, using the addition formula for cos(2(y — x) + 2x) and 
the fact that R(x — y) is even. This is the only essential 
mathematical assumption made on the interaction poten- 
tial R(x — y) in obtaining this family of exact solutions, 
including the normalization condition. The solutions to 
GP theory (|^) found in Il3] are easily recovered by letting 
R(z) = 5(z), i.e., = 1. 

From Eq. (|Tl|d) and P — Vq/ aA(k) , it follows that mea- 
suring the amplitude A(k) of the oscillations for varying 
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external potential wavelengths k leads to the construc- 
tion of the Fourier transform of the interaction potential 
R(2k). Thus by inversion 
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cos(2fcx) 



dk. 
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In principle, this gives a method to determine the pair- 
wise particle interaction potential experimentally. 

Now, we investigate numerically the stability of solu- 
tions of both the local (||) and nonlocal (Q) equations 
without applying any additional perturbations. The 
computational procedure used is a 4th order Runge- 
Kutta method in time and filtered pseudo-spectral 
method in space. Our objective is to explore the question 
of asymptotic equivalence of stability (AES): as the non- 
local equation (Q) approaches the local equation (g) does 
the dynamics of the solution (^) of the nonlocal equation 
converge to the dynamics of the solution (||, with (3 = 1) 
of the local equation? 

For a sufficiently high offset value B, the local solu- 
tion is stable 13 1^] as illustrated in Fig. [j]. In this case, 
the time of onset of instability t* — oo, even under the 
influence of large perturbations Jll|[lJ|] ■ 

To investigate the nonlocal behavior, a choice must be 
made for the interaction potential R(z). A reasonable 
first choice is a Gaussian profile 



R(z) = 



1 



2tt€ 



(13) 



where as e — > 0, R{z) —* 8{z) and the solution (|^) ap- 
proaches the local limit. Thus examining AES amounts 
to examining t* as a function of e. The unstable dynam- 
ics and its spectral evolution are illustrated in Fig. ^ for 
e = 0.01. For e S [0.0025 0.16], numerically t* « 10.1, 
independent of e (see Fig. [|). This result contradicts the 
expectation that t* — > oo as e — > and suggests the pres- 
ence of a beyond-all-orders phenomenon. The Fourier 
spectrum in Fig. provides the primary diagnostic for 
studying the instability and its convergence. Not only 
does it provide the value of t* , it also illustrates that the 
spectral bandwidth, which is the support of the unstable 
modes, is independent of e for e S [0.0025,0.16]. Hence, 
AES is not obtained, even in a convergence-in-measure 
sense. These results raise the following questions: is lack 
of AES due to the choice of R(z)l Is AES a consequence 
of nonlocality or can it occur for local generalizations of 
NLS? 

The choice of R(z) is addressed by choosing differ- 
ent forms of the interaction potential, which are even 
and decaying at infinity. Three additional R(z) choices 
were considered: R(z) = e/(ir(z 2 + e 2 )) with t* = 10.2, 
R(z) — l/(2e) for x <G [— e, e] and R(z) = otherwise with 
t* = 10.1, and R(z) = l/(2e) exp(-|z|/e) with t* = 10.4. 
Thus the lack of AES appears to be a universal feature 
that is independent of the interaction potential. 




FIG. 1. Stable evolution of the stationary solution (|9|) 
with a = 1 (repulsive), (3 — 1 (local), k — 1, Vq = —1, and 
B = 1. 
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FIG. 2. Unstable evolution of the stationary solution ^ 
with a = 1 (repulsive) , k = 1, Vo = —1, B — 1, and 
e = 0.01 in Eq. (fL3h . The bottom figure is a density plot of 
the evolution of the Fourier modes. 



A local generalization of the NLS equation (0) is 

^40 + ^( lrf+e «W (14) 
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where e = z 2 R(z)dz. This equation is derived by a 
change of variables y — x — z followed by a Taylor ex- 
pansion of | ip\ 2 about z = in (]?]). Finally assuming 
R(z) — S(z) results in (pT|). Such an equation has been 
considered in the multidimensional, attractive (a = —1) 
case as a successful means to arrest collapse and blow-up 
of solutions Jl5|,[l6| . The effect of this additional term on 
the stability of solutions was never addressed. A fam- 
ily of exact solutions for this case is given by ([)]) with 
(3 = 1 — 4efc 2 . The dynamics is found to be unstable 
as illustrated in Fig. ^ for e = 0.02. More importantly, 
lim e ^o t* = oo as shown in Fig. |^. In contrast to the 
nonlocal Eq. (0), the presence of the local perturbation 
in Eq. ( |l4|) does not destroy AES. This suggests that non- 
locality is responsible for the beyond-all-orders failure of 
AES. 

Nonlocality for NLS-type equations is a generic feature 
arising in physical systems. Although nonlocality has 
previously been used to prevent the non-physical features 
of collapse and blow-up, its effect on the stability of so- 
lutions appear detrimental. We have demonstrated this 
beyond-all-orders singularity arising from nonlocal per- 
turbations. The specific model considered is the Gross- 
Pitaevskii equation with a standing-light wave potential 
for which an exact family of stationary solutions is con- 
structed. Asymptotic equivalence of stability (AES) for 
the nonlocal equation (^) is not achieved since the nonlo- 
cal equation dynamics does not approach the local equa- 
tion dynamics as e — > 0. This is a truly nonlocal and 
universal phenomenon as illustrated by the study of a 
local correction model and the consideration of different 
interaction kernels. The instability discussed is similar 
to a parametric instability. It differs from it in that it 
is driven by the change of a function, not one parame- 
ter. Further, the equations that are being compared are 
asymptotically equivalent, as are their solutions. There 
is a caveat to our conclusions: in a physical setting, many 
additional effects are present which are excluded by the 
NLS-type model. The success of such models can be ac- 
counted for by the presence of these effects which may 
counteract the instability mechanism found here. 
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FIG. 3. Unstable evolution of the stationary solution ^ 
with a = 1 (repulsive) , k = 1, Vo = — 1, 5 = 1, and 
e = 0.02 in Eq. @. 
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FIG. 4. Time of onset of instability t* as a function of e. 
Note that — log 10 e gives the position of the most significant 
digit in e. 
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